Monthly spatial dynamics of the Bay of Biscay hake-sole-Norway lobster fishery: an ISIS-Fish database

We propose a database to describe the Bay of Biscay mixed demersal European fishery over the period 2010–2020 for the ISIS-Fish simulation tool. It was built upon national and European fishing databases, scientific survey datasets, models estimates, literature, and the formulation of assumptions. It accounts explicitly for spatial and seasonal processes, and for mixed fisheries phenomenons. This database encompasses population dynamics for 3 stocks, hake, sole and Norway lobster, exploitation dynamics for numerous fleets and métiers, and the description of current fishing management, as well as fishers adaptation. A calibration procedure was designed to ensure consistency between all these diverse and heterogeneous parameters compiled and estimated in the ISIS-Fish database and to ensure the model reproduces closely historical catch patterns. This database is a starting point towards operational simulations, of use for understanding the functioning of the fishery, the assessment of management strategies, or delivering short and long-term scenarios. It can be used to reproduce historical catch patterns, with room for improvement on some inter-annual and spatial dynamics.


Detailed model description
The model parametrisation is based on the one used for the projects Myfish (Worsøe Clausen et al., 2016), Benthis (Sala et al., 2014), and Coselmar (Pardo et al., 2017). It was altered to fit our modelling purposes on the hake-sole-Norway lobster mixed fisheries.

Time periods
Selected time period is 2010-2020. The base simulation was calibrated on 2010-2012 (see details section 2.5), validated on 2013-2016, the period of interest to assess current management (including landings obligation) being 2016-2020.
• Sole: only 1 zone (Vigier (2018)'s a priori zone, see figure S1e), given spatial artifacts were appearing in the modelling with 3 zones (all individuals ended up concentrated in only 1 zone; the complexity brought by the 3 zones did not allow us to get more knowledge on the fishery) Seasons Same seasons as Vigier (2018), basing on maturity, movements and growth. For each stock*season*group combination, an accessibility parameter is defined and estimated in the calibration procedure (see section 2.5).
• Hake: 4 seasons/trimesters, at the beginning of which the following events occur. First trimester = aggregation of mature individuals on the shelf break ; Second trimester = Mature individuals begin to spread over the shelf; Third trimester = Mature individuals still present on the shelf break finish to spread over the shelf; Fourth trimester = migration from northern area to Celtic Sea (see areas paragraph "Zones"), also no hake recruitment during that season • Sole: no seasons • Norway lobster: 8 seasons are defined; 1/January, begins with the annual recruitment. Females are inside their burrows, less accessible; 2/ February-March females are inside their burrows, less accessible; 3/ April: Spring moulting, females are more accessible; 4/ May-June, females are more accessible; 5/ July-August, females are more accessible; 6/September, females are inside their burrows, less accessible; 7/ October: Autumn moulting only for immature females and all males, females are inside their burrows, less accessible; 8/ November-December, females are inside their burrows, less accessible (Conan, 1975) Natural mortality • Hake: constant, M=0.5, instead of the more commonly used M=0.4 (ICES, 2014a). We choosed M=0.5 given that the biomass, and then catch in Bay of Biscay had a too string increase across time steps with a lower M. A higher M value would have worsened too much the abundance collapse pattern in Celtic sea and northern area. M=0.5 also proved to provide the best fit to hake catch, without trtiggering an unrealistic discarding pattern.     Table S3: Hake global recruitment ans its spatial and seasonal allocation on 2010-2020 • Norway lobster: M=0.3 for all males and mature females (strictly over 22mm), M=0.2 for all immature females (under or equals 22mm) (Morizur, 1982) Weight-at-class • Hake: W = a * l b , with W the weight in kg, l the mid of a size bin in cm, a = 5.13 * 10 −6 and b = 3.074 (ICES, 1991) • Sole: see table S1 (ICES, 2013) • Norway lobster: W = a * l b , with W the weight in kg, l the mid of a size bin in cm. Males: a = 0.39 and b = 3.18. Females: a = 0.81 and b = 2.97 (Conan, 1978) Maturity • Hake: Maturity ogive from ICES (2010) M l = 1 1+e −0.2(l−42.85) , with M l the proportion of mature individuals inside a size class, and l the mid-size class in cm.
• Hake: at the beginning of each month in January-September. Recruitment is allocated in each zone (Bay of Biscay recruitment zone, Celtic Sea and northern area) and month basing on ICES (2017) and Vigier et al. (2018) estimates. It is then deterministic and independent from biomass variations in ISIS-Fish. Individuals recruited at the beginning of year y, month m and zone z R y,m,z are derived with the following formula: R y,m,z = 1 3 R y * p s,y * p z,y   (2017) with R y the global annual recruitment for year y, p s,y the proportion allocated to the season s (to which the month m belongs) at year y, and p z,y the proportion allocated to the zone z. All values are summarized in table S3. • Sole: at the beginning of each year. Individuals recruit at age 2, ages 0 and 1 not being modelled. On 2010-2016, recruitment is forced to equal ICES (2017) estimates (see table  S4). A Hockey Stock relationship is then used for 2017-2020: 1.8247 * 10 7 , SSB ≥ 9596t (2) • Norway lobster: recruitment is annual, modelled with a Beverton-Holt relationship: R y = αSSBy β+SSBy , with SSB y the mature female SSB in kg at the beginning of year y, R y the recruitment in numbers, α = 1220695327.3 and β = 3866855.3 (ICES, 2014b)

Growth
• Hake: growth increments. Growth modelling is based on ICES (2010)'s assumptions, and is similar to Drouineau (2008);Vigier (2018), except for differences with Vigier et al. (2018) to account for the different widths of size bins: -A linear function between ages 0 and 0.75 with t the age in months,L 0.75 the mean length at age 0.75 year.
von Bertalanffy (von Bertalanffy, 1938) from age 0.75 : with L ∞ the asymptotic length , K a growth rate, and t 0 an "artifact" coefficient supposed to represent the age in months when mean size equals 0. t 0 is the solution of L 0.75 = L ∞ (1 − e −K(x−0.75) ). Following coefficients were estimated by ICES (2010) : L ∞ = 130 cm, K = 0.177319 an −1 ,L 0.75 = 15.8392 cm. From these comes that t 0 0.01727 an.
ISIS-Fish structure being in length only, it is not possible to use directly that relationship, function of the age. Growth is then modelled with: with N t the vector if abundance for all length classes at time step t, and G the transition matrix indicating the probability to go from a length class to another one. G is constant across time steps. Growth increments approach was used to derive G. A mean growth increment for a fish of size l ∆ l (l) in time interval ∆ t is linearly defined beforeL 0.75 . After that length, it is defined with Fabens (1965)'s von Bertalanffy reformulation: Interindividual growth variability is modelled with Drouineau (2008)'s assumptions. Growth increments during a time step are X random variables, of which the mean is derived using von Bertalanffy equation, and of which variance scales to the square of the mean. Then: with α a variance constant (strictly superior to 0). We assume that CV = 1, which gives α = 1.
Finally, it is assumed that growth increments are gamma distributed: X ∼ G(α, β l ), of which density function is : The transition matrix G is then derived as follows: with L the set of length classes, g temp ij the probability to go from length class i to class j during time interval ∆ t , L j and U j respectively lower and upper bounds of length class j, L i and U i those of class i and m i = L i + U i 2 .
However, given that length classes widths were different, the derived probabilities lead to a mean increment being too high just before a change of length class width. g temp are then "corrected", by increasing the probability of not growing and decreasing the probability of growing: with X ik the actual increment when going from class i to class k.
G is then : with n the number of length bins.
• Sole: individuals get older • Norway lobster: growth is modelled by growth increments, similarly to hake, except g ij = g temp ij , since the width of length classes is constant for Norway lobster. One matrix is derived per sex, and to mime the moulting phenomenon, growth is possible only in April (all Norway lobster) and October (males and immature females only). Parameters are : L ∞ = 76mm, K = 0.14year −1 (males and immature females); L ∞ = 56mm, K = 0.11year −1 (mature females) (Conan, 1978) Migrations Defined only for hake, and summarized in table S11. 2 sets of migrations are defined: • One related to spawning and recruitment in the Bay of Biscay: aggregation at the beginning of the year on the shelf break to spawn, and then dispersion on the shelf (Casey and Pereiro, 1995;Guichet, 1996;Poulard, 2001;Alvarez et al., 2004;Woillez et al., 2007). Aggregation occurs at the beginning of January only for mature individuals (except first time step, as mature individuals are already in the spawning zone), and dispersion occurs progressively in April, then July (all individuals that did not spread yet spread). Also, from age 1 (around 20cm), individuals in recruitment zone spread in interim recruitment zone, to model a diffusion towards areas neighbouring the nursery area, at the beginning of each time step.
• One related to an inter-box movement, from the northern area to the Celtic Sea, estimated by Vigier et al. (2018) Accessibility One accessibility parameter per stock*group*season is defined, and estimates with the estimation procedure described later (see section 2.5).
Abundance at initial step For each stock, abundance is provided for the beginning of January 2010 for each zone and class prior to the simulation:   (18) and (19) for functions formulae -Spawning area: all matures individuals over 20cm are assumed to be in that area. For each bin, the number of mature individuals is derived with the maturity ogive (see above) -Recruitment area; all individuals of which length is strictly under 20 cm are assumed to be less then 1 year old and to be recruits, so they are all in this zone -Interim recruitment area: all individuals that have not been allocated to one of the 2 areas above -Presence area: empty Finally, since growth is not modelled for the first time step in ISIS-Fish, the transition matrix has been applied to all zones for hake.

Exploitation model
Gears Gears, their set of species, and their standardisation coefficients are introduced in table  S5. Coefficients values were derived in Marchal (2005)

Selectivities
• Hake: all hake selectivities were derived by Vigier et al. (2018). See table S6 for types of formula and parameter values for each gear. 2 types of formula were used (Methot and Wetzel, 2013), see equations (18) and (19).
• Sole: trawls, Norway lobster trawls and nets catch sole with a constant selectivity of 1. Hence, accessibility is mix of accessibility (catchability sensu Seber (1982)) and selectivity in the case of sole  • Hake: all hake retentions were derived by Vigier et al. (2018). See table S7 for types of formula and parameter values for each gear (only for discarding gears). On top of that, all fleets, including the catch-forced ones, discard their undersized individuals (hake size < 27cm) • Sole is assumed not to be discarded, except because of management constraints. Hence, retention is constant and equals 1 for sole (undersized soles are not modelled) • Norway lobster is assumed not to be discarded, except because of management constraints. All under-sized Norway lobsters are discarded (size < 26mm CL)

Métiers
Here, a métier is the combination of a gear, a zone and a mix of species. The set of modelled métiers was defined: •   Other fleets Spanish longliners and gillnetters and non Bay of Biscay métiers could not be explicitly described by effort due to lack of data. They are instead modelled by a catch in numbers, derived prior to the simulation, to be removed from abundance (in the limit of available individuals) in the post-action of each time step. They catch only hake. See Vigier (2018)'s appendix C.2 for a comprehensive description of these métiers .
That type of modelling for these métiers , combined to the inability of a size-structured model to follow hake cohorts, is thought to be what causes the abundance collapse pattern in the Celtic Sea and the northern area.
Target factors Target factors are multiplicative coefficients in Pelletier et al. (2009)'s eq. (21). In this parametrisation, target factors are the products of 2 components: • A fixed component, for each species*métier *season, derived during SACROIS French effort data analysis (which is the product of a métier component and a season component  MSY transition scenarios from 2016 to 2020 were considered, but had to be dropped out because of our unability to properly estimate F M SY reference point for hake (collapse pattern out of the Bay of Biscay).
Minimum Landing Size All catch has to be discarded (unless landing obligation is constraining) when individuals are smaller than the following lengths: • Hake: 27cm • Sole: undersized individuals are not modelled (under age 2) • Norway lobster: 26 mm CL Landings obligation Landings obligation is implemented form 2016 for hake and sole, with de minimis exemptions for some métiers catching them: • Hake trawlers only: 7% in 2016-2017, 6% in 2018-2019, 5% in 2020 • Hake longliners and gillnetters: no exemption • Sole trawl métiers : 5% ; net métiers : 3% Fishers behaviour Fishers behaviour is modelled by 3 main sets of assumptions: 1/ their reaction to the mix of TAC/TAL, minimum landing size and landing obligation management 2/ Adaptation of fishing effort when TAC/TAL is nearly reached 3/ Adaptation of fishing effort to the current TAL/TAC to spread the catches over the whole year.
1/ At the beginning of a time step, a métier is allocated to a category according to the decision tree figure S2. That tree is explored for each of the 3 stocks, and the highest category number among the species caught by the métier is allocated to the métier . Depending on the category (and "must land" flag for hake), the métier : • Category 0: operates without constraint Figure S2: Decision tree to allocate métiers into categories. The tree is explored by each species*métier combination; each métier is then allocated to the highest category the species it can catch reached • Category 1: can only discard its catch for species whose TAL/TAC is reached; no constraint on other species • Category 2: is forbidden, and its effort is re-distributed to other métiers (see below for details on re-allocation) • Hake case: if the TAC is not reached by the landings obligation is implemented, the métier has to land its catch if the de minimis exemption is reached, operate without constraint otherwise Category 2 métiers effort is then re-allocated in the following order, in each strategy: • to non-forbidden métiers among the same strategy and using the same gear • if said métiers do not exist, to other métiers in the same strategy having a month-strategymétier proportion (see table S21) strictly over 0 • if said métiers do not exist, to other métiers in the same strategy having a month-strategymétier proportion (see table S21) equal to 0 • if said métiers do not exist, the effort is not re-allocated In all cases, effort for the forbidden métier is set to 0.
2/ To avoid an unrealistic discarding peak the month the TAC is reached, the effort is altered when total landings/catches are close to the TAL/TAC value (TAL: hake and sole before 2016, Norway lobster; TAC: hake and sole from 2016). Details are provided in appendix A.3.
3/ To avoid an unrealistic early fishery closure, the catch is spread all over the whole year, using an expected catch for the year to come. Details are provided in appendix A.3. 2/ and 3/ are based on nominal effort times the inter-annual variation factor for the time step. Hence target factors per se are not included in these calculations.

Calibration
Accessibilities and targeting factors were calibrated on catch observations and fishing mortalities estimates, following a procedure exposed below.

Hake accessibility
, with F O 1 , F O 2 and F O 3 the 3 objective function components (first component focusing on length composition for each season and fleets ; the second focusing on catch in weight for each season and fleet; the third focusing on annual catch in weight); ω LF D .,.,.
weightings on weight observations. See Vigier et al. (2018);Vigier (2018) for the rationale behind the weightings derivation; C sim .,.,.,. catch simulated in ISIS-Fish; C obs .,.,.,. the observed catch. α, β et γ weightings on objective function components, respectively10, 1 and 35 (empirical), allowing the balance the influence of each objective function component, and favouring the third component, to improve the ability to reprocude annual catch in weight.
The calibration procedure consists in computing accessibility values minimizing the objective funtion using an iterative algorithm. At the beginning of the initial iteration, accessibilties values were set as inVigier et al. (2018). A batch of simulations was run, each using a multiplicator value applied to all hake accessibilities, varying between 0.5 and 2.5, with a 0.05 step. The objective function was derived then computed for all this values combinations and the optimal multiplicator value was identified for each parameter (quarter).
A new accessibility value was deduced, Next iteration is performed by multiplying the accessibility value of each quarter by their respective optimal multiplicator value.
Next iterations were initiated similarly as what was described in last paragraph. We stopped the estimation procedure when after an iteration, accessibilities values were not changed. Final values are provided in table S17.

Sole accessibility
, with A the set of age classes, Y the set of years, F obs . the "observed" working group estimate, F sim . the ISIS-Fish estimate. At the beginning of the initial iteration, The calibration algorithm starts with accessibilties values set with values in the original database (V0). A single simulation was run, and the objective function was then computed. If it was deemed too high (i.e. F obs a,y F sim a,y too far from 1 for at least one age class), all accessibilties values were multiplitied The following tested values of sole accessibility is derived by multiplying the current accessibility values by F obs a,y F sim a,y ratio. Next iterations were initiated similarly as what was described in last paragraph.We The procedure was stopped the estimation procedure after 3 iterations, as all F obs a,y F sim a,y ratios were close to 1 (subjective). Final values are provided in table S18.

Norway lobster accessibility
, with q a quarter, L the set of length classes, S the set of sexes, C obs .,. the observed catch in numbers, C sim .,.
the simulated catch in numbers One estimation procedure was run for each quarter, in the order of quarters, as the catch done in the first quarter can influence the catch done in further quarters.
At the beginning of the initial iteration, accessibilties values were set as in the initial database (V0). A single simulation was run, and the objective function was then derived. If it was deemed to high (i.e. Target factors T arf met,sp,q,y drive how the effort is distributed between population sp, métier met and season*year (q,y) combinations. They were split in 3 components: a fixed component T arf met,sp,q for each métier met, population sp and quarter q, derived from the effort dataset analysis, another fixed component T arf gmet,sp,q,y , for each group of métiers gmet sharing the same gear (longliners, gillnetters, roundfish trawler (coastal), roundfish trawler (not coastal), Norway lobster trawler), population sp, quarter q and year y, driving inter-annual variations of fishing effort, derived from catch observations (no quarter dimension for sole and Norway lobster for that component), and finally an estimated component, T arf est gmet,sp,q , allowing to tune the model's dynamics to observed catch. This section focuses on the estimation of teh latter. The 3 components were combined as follows: T arf met,sp,q,y = T arf met,sp * T arf gmet,sp,q,y * T arf est gmet,sp,q As for accessibilities, target factors estimated components were estimated jointly within a population, but not between populations.
Hake target factors 20 parameters were defined, for each combination of the 5 five group of métiers (longliners, gillnetters, roundfish trawler (coastal), whitefish trawler (not coastal)) and 4 quarters. We used the same data and objective function as for hake's accessibilties estimation. At the beginning of the first iteration, each estimated component was set equal to the ratio of observed catch over simulated catch C obs gmet,sp,q C sim gmet,sp,q . A set of 21 simulations was then launched, each of them with values of estimated components of the target factor varying between 0.5 times and 1.5 times the latter ratio, with a step of 0.05. After the runs, objective functions were derived, and an optimal ratio value (minimizing the objective function) was identified. If theses values were 1 for all target factors, the procedure was stopped. If not, a new iteration of this procedure was done, by defining the estimated components of target factors being equal to Sole target factors : 1 estimated component per group of métiers (gillnetters, Norway lobster trawlers and whitefish trawlers) and quarter. We were provided with sole catch in weight on 2010-2012 for each métier and quarter. Just after hake's target factors estimation, a simulation was run, and for each group of métiers and quarter combination, the average of yearly ratios of observed over simulated catch of sole C obs gmet,sp,q C sim gmet,sp,q was derived. These ratios were used as multiplicative of target factors for sole (T arf est gmet,sp,q ). After a simulation using these latter target factors, all the latter ratios were close to 1, showing the catch distribution between métiers and quarters was more realistic. Hence no further iteration has been done. Estimated values are provided in table S29.
Norway lobster target factors : 1 estimated component per group of métiers (Norway lobster trawlers and whitefish trawlers). We were provided with monthly Norway lobster landings data per length and sex class for 2010. Just after sole's target factors estimation, a simulation was run, and for each group of métiers, the average of monthly ratios of observed over simulated catch of Norway lobster C obs gmet,sp,m C sim gmet,sp,m was derived. The average of these ratios for each group of métiers were used as multiplicative of target factors for Norway lobster (T arf est gmet,sp,q ). After a simulation using these latter target factors, all the latter ratios were close to 1, showing the landings distribution between métiers was more realistic. Hence no further iteration has been done. Estimated values are provided in table S30.
A Additional information about ISIS-Fish parametrisation

A.3 Management
Reducing effort and discards when the TAC is nearly reached Condition to consider if TAC is near-reached during year y at the beginning of month m (numbering beginning at 1): the expected catch this month (the monthly mean over the current year) is bigger than the remaining catch.
The above condition is verified each time-step from 2013, before métier flagging, until the TAC is considered near reached or reached. The catch is either the landings or the catch, depending on whether the TAC is a TAL or a TAC (always a TAL for nephrops; TAL for hake and sole until 2015, TAC then). If the TAC/TAL is considered near reached at step m, it is automatically considered reached at step m + 1, whether it is actually reached.
Ratio to apply to the month-strategy-métier proportion: One ratio is derived per species whose TAC is considered near reached. For each métier catching several of these species, the samllest ratio is to be applied to the corresponding strategy-métiermonth proportion.
Once the TAC is considered reached, only if the landing obligation applies, a few more calculations are done in order to limit the discards during the first month of reached TAC (reducing an artifact). these are done just after the re-allocation of forbidden métiers effort.
Computing the allowed discards this month: Computing the ratio to be applied to an effort value: are not derived on the same duration, because almost no catch happens after the TAC is reached and/or near reached, and also beacause using m instead of m T increases r disc , and hence worsens the discarding artifact. This ratio is computed for each species. Métiers catching several species use the smallest ratio. A maximal strategy-métier -month proportion is derived as following: StrM etM onth max str,met,m = r disc,sp m T t=1 StrM etM onth str,met,t The strategy-métier -month proportions are then capped: This capping occurs after effort has potentially been re-allocated to métier met.
Ideally, this capping should ensure that the de minimis exemptions are respected, but given q, T arf and biomass variations, more discards occur. It is however still less than before these calculations. Also, some métiers may have no effort allocated to them at the beginning of year; but not at the end. Hence, they are forbidden, and maybe re-allocated effort late in the year, explaining why marginal discarding may occur months after the TAC was reached.
Adapting effort to the TAC value Since the previous correction is not enough to eliminate all the artifacts, the following has been implemented. The idea is that the effort is altered for the whole year once the TAC has been derived to prevent extreme overshooting or undershooting (done once per year, in January, from 2013). Here catch and TAC are used irregardless of working on landings and TALs (hake and sole before 2016, Nephrops) or catch and TACs (hake and sole from 2016).
The idea is to compute for each year for CPUE in the Bay of Biscay for each stock: CP U E BoB,s,y = C BoB,s,y Ey Effort is only accounted for métiers catching the stock s. The catch of Spanish longliners and gillnetters is accounted for hake, but not their effort. Hence, for hake, it is a catch per unit of effort not coming from the Spanish longliners and gillnetters; still, this is enough to correct the effort.
For each stock, a mean CPUE over years is then derived: CP U E BoB,s = y CP U E BoB,s,y num years (27) From which an expected catch is derived, using what we already know on expected effort for the current year (equals 2010 effort, times the inter-annual variations of effort for the current year, at the gear level): C exp s,y =CP U E BoB,s * E(y) Then a ratio is derived: And the inter-annual variation components of the target factors are accordingly altered: V arsEf f N EW met,s,y = r s * V arsEf f OLD met,s,y For each métier met, the ratio r s is: • The smallest of all ratios if at least one stock is expected to be overshot (at least one ratio is smaller than 1) • The biggest of all ratios if none stock is expected to be overshot (all ratios bigger or equal than/to 1)